DOCUMENT RESUME 



ED 218. 126 

AUTHOR ' 
TITLE 



i 



SE 038 236 



INSTITUTION ' 
SPONS AGENCY 
PUB DATE 
GiiANT 
NOTE 

EDRS^ PRICE 
DESCRIPTORS 



Rheinboldt, Werner c'. ' - ^ - 

Homer • s .Scheme and Related* Algorithms. Applications 
of Computer Science, [and] Algorithms for Finding 
Zeros of Functions. Computer Science. Modules and 
Monographs in Undergraduate Mathematics and Its 
Applications Pr^jkct. UMAP Units 263 and 26^4. 
EflucatFon Development Center, Inc., Newton Mass . 
National Science Foundation, Washii^gton, ' D.C. 



National 
80 

SED-76-196i5-A02 
59p. 



MFOl Plus Postage. PC Not Available \frOm EDRS. 
♦Algorithms;. *College. MatHematics; Computer Programs; 
♦Computer Science; *Computer Science Education; 
Higher Education; Instructional Materials; *L'eQrning 
Modules; *Mathemat ical Applications;. Problem Solving; 
. Supplementary Reading Materials 
IDENTIJ'IERS *BASIC Programing Language 

ABSTRACT , , ' ' 

This material contains two units which view, 
applications of computer science^ The first of these exajnines 
Horner's scheme; a^io is .designed, to injstruct the user on how io apply 
both this scheme and' related algorithms. ,The second ufnit aims for 
student underjstanding of standard bisection, secant, and Newton 
methods of root finding and appreciation of their limitations and 
strong ppints.tAn i^ntroducti*n to more recent root finding methods is 
also provided^ Both modul^ contain exercises, and answer^s to these 
problems are given at^he conclusion of ^ each uni^t. (MP) ^ 



-r 



****** Y**** ********* ************* 

* Reproductions supplied by EDRS are the best that Ulin be made 

* irom the' original documientl 

********* ^***^******* **************** *************'^^^^^^^^^^^^^^^^^^^ 



oo 

Q 
UJ 



vS 



UMAP 

MODULES AND 
MONOGRAPHS IN 
UNDERGRADUATE 
MATHEMATICS ' 
AND ITS . 
APPLICATIONS 



^ ^ ^ 



CS3 



CS! 



X 

>^ O >^ O X <i> 
o > 



^ :r 



X 

t 

Q 
Co 



> o 



lU 



X 
hi 



MODULE 



263 



Horner's 
Scheme and 
Related 
Algorithms 

by Werner C. Rheinboldt 















1, -2' 


p 


1 


-1 




0 . -1 


3 




2 


q,{x) * 


1 . -3 

1 

'O -1 


3 
4 


-2 
-7 


1 =P( 


•q,(x) 


1 -4 


7 


-9 






0 -1 


S 






q^Cxi 


1 ^-^ 
0 -1 


12 


1 

" 2! 


p;'(-i) 


qjx) 


1 1 -6 


' 1 
" 3! 


p- (- 


1) 











Applications of Computer Science 



p • ? 



EOUC^I 



U.S. DEPARTftllENT OF EDUCATION 
M NATIONAL INSTITUTE OF EOUCATION 

EOUCi^TIONAL RESOURCES INFORMATION ^ 

CENTER (ERIC) . " 
document has been reproduced as* 
BirkhdUS6r Boston Inc. fweived from the person or organization 

OCW^ O* ongtr^ating it * • ^ . 

3oO Green street Mmor changes have been made to improve "PERMISSION TO REPRODUCE THiS 

' Cambridge MA 02139 reptoducvor^ quauw material in microfiche only 

' ' HAS BEEN GRANTED BY 



• Points of view or optnions stated in this docu 



ment do not necessarily represent official NIE -/^ - 0 , . 

'pos.tK)n or policy I UiAJju)n(U ^>t3^/tCg^ 



TO THE EdJCATIONAL RESCJUpCES 
INFORMATION CENTER (ERIC)." 



intevmodulca* DeecHption Sheet : UMAP Unit 263 

' Title : HORNER'S SCHEME ANO RaATED ALGORITHMS 

^uthor: Werner C. Rhelnboldt 

Department of Mathematics 

and Statistics 
University of Pittsburgh 
Pittsburgh,^ 1526] 



Review Stage/Date :'* IV n 7/30/80 

Classification : ' APPi, COMP St I ^ 

Prerequisite' Skiyls : ^ 

Def i nition of a derivative. 

2. High school algebra. • ' 

3. Fundamentals of how t<? draw a flowchart. . 

Output Ski'lls : - ^ — ^ 

1. Be able to apply Ho^'ner's a^d the related algorithms of the 

The Project would like to fhank C^arlBs Votaw of Fort Hays 
University arvcj Brian' J. Winkel of Alb^ion dollege for their Vevi 
and all others who assisteti in the production of this unit. 



uni t . 



State* 

ews. 



^ ^ This material was fiel-d- tested and/or student reviewed in, 
Jiminiry f6rm 'i n class ^b^ J<5natfian ChbfftaJ'of The Grot|R School; 
,B. .M;hase of ^ssi^h to4.1ege T Dfona Id G. Halm oT Oalcland Universi 
Mary Haywand'6f Rbbef ts^ Wpsteyan College, and has been revised 
bas'is of d§ta. recei^ve4 frdm» these s-ftes^. 

^ ' ' . - . \ 

Tbt s.matferial was prepared with tbe partial support of Nat 
Stt^nce Foundation Grant* No. ^ft76- 1^6 15 Aa2. Reconmendatjons 
*'^pr*essed are th6se.;of the author and do not necessarily rfefVect' 
> -views^;qf ^th^ NSf. or ithe'copyrighl holder. 



.pce-'. 
^Gerfe 
ty,^n^; 
on Bhe 



6r1al 
exT . 
tfie 



3 • 



©J980 ^DC/Protect" UMA.P 
Al 1 rights jresej'ved. 



ERIC 



HORNER'S SCHEME AND RELATED ALGORITHMS 

' Werner C. Rheinboldt. ^ \ 
Departme^,nt of* Mathematics - 
and Statistics 
• University of Pittsburgh 
Pittsburgh,' PA 15261 

TABLE OF CONTENTS ' / 

1. INTROOUCTIOK ... ' ' , 

' ■ ' * ^ .... V 

^ . ^ ' ' ' ' ■ ^^-^ . 

2. HORNER'S SCHEME / . ^ & ^ o 

* ^ ' 

3. ^ IMPLEMENTATION OF HORNER'^S SCHEME V . k * * 

CONVEJfelOtf TO DECIMAL REPRESENTATION .... ♦ " ' .5" 

■\'- ■ " . ■■■ 

5. -\ HORNER'S SCHEME AND POLYNOMIAL DIVISIOH 6 

6. ..HORNER'S SCHEME AND THE DERIVATIVES , . . 10 

* ' / " ' ' ' 

^ * ** . • - ► - 

7. OUTLOOK .... • - 

< / 

\ " 

8. ANSWfRS TO EXERCISES ^5 

............. . 

/ ■ ■ 



1> INTRODUCTION ' . - . 

V 

Our basig problem is the computational evaluation 
of a polynomial 

(!)■ P(x) = a„x" . a„.,x-l\ a^.^x"-^ . 
... + a x2 +^ a ,x ? a, 

and its derivatives p ' (x) , p".(x) , p^"^(x) for a 

prescribed. v^alue x = x^,. Here the coefficients a , 
^n*l' ^0 given real numbers. 

/ Let us consider, first a simple cubic polynomial^ 

(2) P(^x) = 3x^ - 4x' + 2x - 3. 

For any given number the ev'aluation of pCxJ does 
not present any principal difficulties. We may compute 
and x^ aiid they combinefthem together appropriately. 
>Ii) an informal prog^ramming language this may be written 
as- thre' following" algOYit;hmi^ ^ ^ . 

1 Input {x^i . \ 

2. u 

(3) 3. ' V 

4 . p := 3v - 4u + ,2Xg - 3 

5. Output: *{Xo,p}. 

, (:= is used to represent assignment of a value to a 
variable.) 

Altogether there are five multiplications and three 
additions (or subtractions). For the general polynomial 
(l)'this approach would require the computations: 

(4) ' r • 

P = W * ^n-'l"n-l * ••• * ^ ^0- ' . 

Thus altogether we need 2n-l multiplications and n 
additions. Suppose that a particular computer uses a 
sec and u sec for any addition (or , subtraction) and 

1 




multiplication, respectively, Th^n our method (4) t'akes 
at* least (2n-l)vi + na seconds, W-ithout question a prac; 
t^ical computer program woulcL run longer than that, sii)ce 
it' takes time to fetch and store the data, to control 
the loop involved in (4), and to-perform the input and 
output. But the overall time*$houlJ be* proportional to 
((2n-l)^i + na) . The next section shows that can*do ,^ 
better than th^%, 

2. HORNER'S SCHEME * 

How can we reduce the numbtfV of arithmetic operations 

» 

in the evaluation of a polynomial*? The clue is a suitable 
^factoring of p(x) . In fact, (2) can be written as follows 

p(x) = ((3x-4)x+2)x-3: 

Now ther|. are only three multiplications and three add'i- 
tion^. That does not appear to be much, of a savings in 
this case but it does represent a big savings when the 
•degree of p(x) goes up. . , 

We sh^ll discuss this approach of evaluatin^a poly- 
nomia],-*in the form of a scheme foi^ hand computatioi>. Let 
a general cubic polynomial be giv^en: ^ » * 

(5) p(x) = ax^ + ax^+ax + a. 

We may rewrite it as _^ ^ 

p(x) =^((a x + a )x+a )x"i-a . 

3 2 10 

To evaluate fhis for x = Xq we use a t^ble with three 
roAvs and four columns. Into^the first row we write the 
four coefficients a^., , a^, a^ (in that order), and 
into column 3, row 2 we place a zerp. The rest of the 
boxes are initially blank. 



, Col 3 'Col%' 


Col 1 


Col 0 


. ' Bow 1 'a(^ 


• 






^ ' (6)' Row 2 • 0 








Row f 83 , 









Tfee computation proceeds from**the left-most cofumn to the 
Vight^an^ corfsists of looping through the folloiving twQ j> 
•steps'^for k = 3,2,1,0:* ^ . , 

• » < 
(a) Add tjie numbers m- rows- 1» and 2 of column k 

and write the result into row 3 of the same 
- column . 

^ (b) For k > 1 multiply the number in row 3 of 
column- k by and pUce the result into 
row 2 of column k-i . \ \ . 

This process is indicatedt by arrows in ^the^able (6) 
and the results are indicated in each fieU. The final 
result in row S^of.the last column is the value of ,p. 
ar^he point x^." This method of computing the value 
of a polynomial is ^called Horner's Scheme.' 

^ We give Horner » s -scheme^f or our ^spetial polynomial 
(4) and two different values of x : 
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(8) 



another example, consider the quartic polyndnli^l 



P(x) = -x" " 2x' + X - 1. 



Note. here that the co^efficient of x^ is zeror it should 
be included, in the computation with *at value. 
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Evaluate (2) at = 1 , = 


\ 

■^ ^0 


= 10, Check your answers. \ 
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^ 2. 


Evaluate (8)' at x^ » 0, x^ =• 




♦ 

« 2, What can you say t 




0* . 


about the behavior of p(x) in 
1 < X < 27 


the i 


ntervals -1 < x < 0 and ' 





3. Differentiate the polynomial (8) and evaluate the resulting* 
cubi.c polynomial at x = -1 and x = 2, ^ 



3. IMPLEMENTATION OF HORNER'S SCHEME 



How can we program (6) for a general polynomial (1)? 
Assume' that, the coefficients a^, a,, ... a are stored in 
an array of length *n+r. If we are in column k of Table 
(6) and the number in row 3 of that column is stored in 
p, then^the numbers In rows 2 and 3 t)f co-lumn k-1 will ' 
J)'e pxj and px^ + aj^ respectively. (In the left-most 
column the corresponding numbers are, of course, 0 and 
a^^.) Thus we c^n write the overall process in the form 
of a simple loop: 

..Input U^y a^ , . . 



1 .• 

2 . 



^n' 



(10) For k = n-1, n-2, 

3.1 p := pXq.+ aj^ 
4. . Output {x'(^,p}. 



0 do 



V 



^ - Each execution of step 3 involves one multiplication 
and one addition, that is, altogether there are n multi- 
plications and additions each requiring n(vi?*-a) seconds. 



3 



J 



This represents a considerable saving^ over the (2n-l)p+na 
.seconds needed for (4). / * , - 



Exerci ses 

^. |raw a flow chart for the algorithm (10). 

5. If a progranxnable calculator or computer ^is availabhe, 
implement (lO) for cubicand quartic polynomials. T^est 
your, program* with the polynomials (3)^ond (8) ^nd the 
values of x used in Exercises 1 and 2 of Section' 2. 

6. Do a hand calculation or run your program on another 
polynomial such as p(x) = x^ + x^ + x + 1 at x =1, 



0 



4. • CONVERSION TO DECIMAL J^EPRESENTATION 

As an application of the Hpr.ner scheme, consider 
the question of finding the decimal representation of 
. some integer N = (a^aj^*. • . . aj^ given nn base b' notation. 
FDr example, let b = 2 and- N the binary number 

(11) .N = 110011. 

Generally, the notation N « Uj^a^^^^ . . .aj^ means that 

Jn other, words, if we. introduce the polynomial (1), then 
we, have N '= p(b) . ^ Thus , in the .case of^ (11); we need to 
evaluate the polynomial ^ ^ 

^ p(x) = .+ x"* + X + 1 . % 

at Xq = 2. The*Horner scheme for this has the form 

X = 2' • 



1 • 


1 
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0" 


. 0 
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• * 1 
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12 


24 


50 
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— w 

» 3 
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12 


25 


51 



and hence our binary number llOOll^has the decimal repre- 
sentaftion 51 . 



Exerci-ses • • 

7. * Convert^the binary integers 10101 ^ 1111, 10000 to decimal 
, representations. ' . ' 

^ 8. Coavert the Integers (7^013)g, ' (l 12ipi)^, (^*^tpl)^ 
• decimal representations.. 



S. HORNEI^^S SCHEME AND ^POLYNOMIAL DIVISION*^ 



We t^etu/n to our ^go.rit'hm (II)), but this 'time we 
retain th& successive values of the variable p in steps 
2 and 3.1, that is, the entri>^$ in ro'w 3 of Table (6). 
rewrite the algorithm as follows: 



1. - Input {a^^, a^,x^} , 

2. p := a^ 



n 



3. 'For k = n-1, n-2-,, . . . , 0 do 
" 5. 1 := p 

. 4. Output {x^,"p, qQ,.--..,iqn.i}K ^ < 

>ience betwee^the coefficients §0, , a^, q^, q^. 

apd p we have jthe relations: . •» 



(13) 



J.- 



'6 
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Evidently in Table (6) , ^ , q^, p are "the ..numbers i^ 
row 3, ' 

< ' Me now introduce the new ^palynomi al * * 

It is related to" p(x) /ria a simple formula', rn'fact, 
using the f6rraulas (li) we f fnd. that . > ' 



^^nd hence bec^ause oV p> pCxJ that 

(15) p(x) =*q(,x) (:j-x^) +.p(x ) '. 

Thus, q(x) is the result of the division of the 
•pol^nomi^r p(xO the linear factor x-^x^ and pCx^) 
is the remainder* Horner»s scheme^is Indeed only a 
slight^ly con^ense^'form.of the standard dgi^fision of a 
polynomial J)y a linear factor! To illustrate this 'we 
writ^ this process in its familiar f6rm for the- cubic 
polynomial (5) (recall = aj): 

x-Xj AjX^ + a^x^ + aix-a^ 

• qz^yqax^x' 

•• • , '■ + ajX 

"qix' -. q^x^x' ^ 
"V * qo^^^o 
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which means that ^ 



V^a^x'+a^x^+a^x+a^ ^ (x -x ^) (q^^x'+q^ x+q^) +p 



in agreement with (15). . . ^ 

r 

/ -More general4-y, this long division may be applied 
^ to divide any polynomial by anonzerb polynomial of 
lower degree. For instance", the division •of (8) -by 
u(x) =*x^-x+2 proceeds as follows: 

^ 

x2-x-3 



x2-x+2 /x'*-2x'+0-x^+x-l * 



9 



(17) - - 2x^,*x ^ ^ 

-x^ ^ x^V2x 
' • " - 3x2'+5x '1 " • ^ 

- - 3x^ + 3x -6 ' ' 



which means that « " ^ . • ^ 

. x'*-2x^ + x-l = (x2-x+2)^(x^-x-33 +5*. ' \* 
- For ;th^ general polynotnial ^1) and any 
(IB) !^(x) = u^x^'+u^^^x"^-^ +u^x+u^/m<n,. u^y(U 

'this' division, algorithnTlnay be;' written in'^the following, 
form: • ' . ♦ 

^ . 1.; input (a,, . a^, t5^, . 

'2. ,For j = p, 1,*..., n do r. := a. 
tf3. For k =^-m, n-m-1^, O de 

„(19) ^ • 3.,i; q, U r^Ju^ " . . 

\ 3.2 FoT. j ='ra+k-l, in+k-2, 

For thfe resultioig polynomiajl^s ♦ . ^ 



•(20) q(x) = qn.;x"«-'^^....-q_^x.q^, r (x-) = '.f^. ^x'"' V . . . *r jx^r „ , 
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we then have i * 

(21) p(x) = q(x)u(.x) +>(x).^ ^ * . 

In the special case of u(c) = x-x^^ that is, m=l, Ui = l/' - 
u^ = -Xjj, the algor,ithm reduces essentially .to (12). The 
only difference is that the remainde^ is now a polynomial 
which we initialize in step 2 as r(x) = p(x). Formerly 
we knew that the remainder is a ^constant which may be 
initialized -a^ p \ a^^. 4 < 

*, 

' The algorithm (19) is 'probably easiest to understand 
by going in detail through the f oil owingt»l!'W»^- hand divi- 
sion process: 

a; 



U, X +u, x+u^ 



^x'^ + rjX^ + rjX^ + r^x+r^ 



r =q u , 



^2^2^ "'^2"l^ *^2"0^ 



r^jX^+r^x^ + r ^x 



q u x^+q u x^+q u X 

^12 ^11 ^10 



r X +r x+r 

2 \ 0 



f r^x+r^ 



remainder 



^Exferctses . * 

9. If^erfbrm the division (16) for the polynomial (2) ^and = 2. 
^..Compar^ your resul rs wi th ,those of (7). " * 

JO. As in (17). divide p(;t) = x^+x^-x'*+29c^-x+2 by utx) ^ 2x^+2x^-x+3. 
Thefv follow the same steps in Algorithm (l9). 

n. If a programmable calculator or cgmputer is^avai l^ble, implement 
(19) for reasonable values of n£m>0. Test your program v^lth.the 
polynomials of,(l7) and Exercise 2 above. 



Mr (Optional) Show that there is only one pair of polynomials q(x), 
r(x) witJ^degree r(x) < degree u(x) that satisfies (21). 



6. HORNER'S SCHEME AND THE DERI VAT I VES 



We return to the basic formula (15) 
tufns out to be the difference quotient 



q(x) 



P(x)-p(x^) 



Since 



'we expect that q(x^) is the val 
P at x„ 



of the derivative of 



we obtain 
(22) 
whence indeed 
(23) ' 



In fact,.;))' applying th^^product rule to (15) 
P/(x) - q'(x){x-x^) + q(x) , 



P'(Xo) = q(xj 



Thus p'(Xo) may be evaluated by applying Horner's scheme 



to q. For our- example (2) and 
follows: 



= 2 this looks as 

s 



2 



3 
0 



-4 

6 



2.^ 



-3 
'12 



(24) 



3 
0 



2 
6 



6 
16 



9 = p(2) 



8 



22 



P^(2) 

To ^implement this, note that each column of (24) can be 
computed from the column on its left-. Thus, we don't 
have to complete the evaluation of p(2),.i.e., fill in 
all of row 3, before finding the valuer of p'(2). How- 
ever, observe also that in the last column^nly p itself 
is evaluated. Thus we mAy extend (10). as follows: 



10 
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1, Input {a,, x ) 



.(25) 



0 » • • • > '*n*' 0 

■ 3. p- := p 

4. For k = n-l,,n-,2, . . . t, 1 d'o 

4.1 p ':= px^o + 

4.2 p'» := p»x^ + p. 

5. p :» pi^ + , ^ . 

6. Output {x , p, p' }. . 

The process may- be exteftded to higher derivatives. 
For this note that the repeated application of the Horner 
scheme results ih*a sequence of divisions: 

='qi(x)(x-x,) + p(xj 

> q^Cx) = q^(x)(x-xj + q^Cx^) 

^26a) ^2^^^ = .q3(x)(x'-x^) + q (X ) 



%M = qj^.i U)(x-x J . qj^CxJ 

where qj(x) denotes ©ur original q(x). At each step 
the degree the q*s decreases exactly by one; that is, 
t.he degree of q^/x) is n-1; for q^Cx) it is n-^ ^nd^. 
generally qj^(x) has degree n-k, Thbs- q^^Cx) is *a'*cons^t ant 
and we have qj^+^Cx) = 0, and our sequence of equations (^6) 
ends with 



• (26b) %.l^^) = qn(^)(^-Xo) ^ qn.l^^o> ^ ^ ' 
q„(x) =-qj,(Xol. ' '\ 

We multiply the k^^ equation by (x-xj^, Ic' = 0, n, 

and add all of them together. Then 'for k = -1, n, the 
^ term qj^(x)(x-xj^ arising on the Jeft of the k^^ equation 
cancels agaii>st the same term on the right in the (k-1)^^ 

equation, Henc^ we obtain ^ 
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x-xj^ ... + 



(27) P(xr^p(xj+q^(x^)(x-xj+q^(xj(x-xj 

9 i 

By differentiating this k tim^s^O < k' < n,' the first 

St * t h 

(k-1) terms disappcfar, the k " term becomes klq^^CxJ^ 
and all subsequent terms still have a nonzero power of 
(x-Xq) as^^factor. Thus for x = x^ these" terms become 
zera-and we find that 



(28) 

Moreover, (27) becomes 



(29) p(_x) = P(xj+p'(x,)(x-xj +^ P*'(xJ(x-xj2+. 

' ^P^"^(Xo)(x^o)". 
This is the Taylor expansion o,f p(x) at x = x^. 

The sequence of divisions (26a/b) is, of course, 
computed by means of repeated applj.cation of the Horner 
scheme. For example, in the case of the quartic poly-, 
nomial (8) we obtain for x^, =^ the following' results: 



-1 



(30) 



P(x) 
qj(x) 

qjx) 

9 

qjx) 



-2 
-1 



. 1 

/Is 



0 



-2 
-^7 



-4 
-1 



-5 
-1 



?(-!) 



P'(-I) 



12" = ^ p"(-l) 



and thus 



12 
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(31) 



p(x) = 1 - 9(x+l) + 12(x+l)2 • 6(x+l)^ + (x+1)' 



Besides fStoviding us with a, simple method for the 

*^valuation of .the 3efivatives of p(x) at a'^given point 
x'^ we -have obtained here also an algorithm for 

rewriting p(x) ip terms of the powers of (x-x^) i^'nst-ead 
of those of x.^ V ^ V ' 

» In extension of 'the algorithm (.25) the entire 
process cai^ be written as follows: 



(32) 



1. Input '{a^ , .. . ^n' ^0^' 

2. For k 1, n do pj^ : = 

3. For k = n- 1,, n* 2 , *. . . , 0 do . 
.3.1 p^^fi.^ + aj^ 

3.2 For j = L, \ , . k do 

' Pj ^= Pj^o t Pj.i. 

4c . Output {x^, p^, PJ'.^^ 



The resulting values are 

and hence are exactly the coefficients of the "shifted" 
polynomial (31) . 

As (2S) the algorithm (32) computes the data column- 
wise from left to^ rrght. The computational process ^s" 
eas.ily followed in |:he next table. • 



Pi 
P3 



a,* ZajX^+Bj * 3ajX^, + 2ajX,+a, 



^3*^ 3ajX^+a, 



a,x,' + ajX^ + a,oc/ a, 



13 
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Exercises . ' 

13. Verify by direct differentiation andUvd-lilat ion of the 
resulting polynomials, the resul tsj g iven in (2^)7 (30) , 
and (31). . . i 

I I'*. Follow in detail,,the algorithm (32) on the example (30). 

15. If a progrlmmable calculator or computer is available, 
implement (32) for reasonable values of n. Test your 
program with the^data in the example {2k) and'(30). 

16. Write out in detail t>he proof of (27), and (28). \ 

17. . Compute the coeff ici.ents of p(x-l^) for 

p(x) = x^-6x^ + ISx"* - lOx^ - 15x2 + i*x - 9 



7 . OUT LOOK 

. , • ^• 

• The basic method jiamed in the'title pf this unit 

'was given by W. G. Horner in the early nineteenth qentury 
liri connection with an. efficient me.thod for finding the 
ciefficients of p-(x-x^) , [Philosophical Transactions, 
Royal Society of London 1*9^819, 308-335]. But'the' 
factorization 

P(x) =•(•••((Vo^^.l)^^an-l)^.^-:0 * , 
on which it is based was already u'S'e'd^y Isaac New^:on 
some hundred years earlier [Analysis per Quant'itatem ^ . 

Series, London, 1711]. - ^ ' " * , 

We saw earlier that Horn^^s method usJs f^wer 
operations than,, f ox 'ins tance , ihe approach indicated 
in (4); It* can be shown that wh^ the inputs to our , 
algorithm (10) are arbitrary constants, that is, when 
we are not using /ny further information about them*, 
then there is no other algorithm which computers p with 
less- than n multiplications a^d n additions. * 



14 
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- In practice the computations in all our- algorithms 
are performed in floating point arithmetic on some 
computer. ' Then round-off errors are introduced and the 
question arises how they affect the results. For in- 
stance, it turns out that with .increasing |x*|, (absolute 
value of x^), the result *of the Horner algorithm (10) 
may be increas^ingly inaccurate. The situation is more 
complex when it comes to the other algorithms given here.* 

8.^ ANSWERS TO EXERCISES 



t. X = 1 

0 



i 4 
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X = -1 
0 « 



X = 10 

0 



f. X «'0 

0 



X = 1 

0 



* 3 




2 -3 


i 

t 


0 


i 


-1 • 1 




3 


-1 


1 -2 


P(l) = -2 


3 


-k 


2 -3 




0 


-3 


7 -9 




3 


-7 


9 -12 


p(-l) = -12 


3 


'k 


2 ^3 




0 


30 


260 2620 




3 ^ 
' 1 


26 

-2 ■ 


262 2617 
0-^ 1 


p(lO) - 2617 

rl 


0 


0 


0 0 


0 


1 


-2 


0 1 






-2 


0 'l 




0 


1 


-1 . -1 


0 


1 


-1 


-1 0 


— v^ • 

-1 ■ ■ _ p(i) = -1 








* 






% 
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X- - 2 
0 > 



2^ 



P(2) ^ 1 




-1 < X < 0 Since p(-l ) > 0, p(0) < 0 
^ p(x) must cross the^x-axls, 

i.e. 'have a zero, between 
, -1 and 0. 

^ 1 1 X < 2 Since p(l) < 0. p(2),>0 

p(x) must cross the'x-axis 
between 1 and 2.* 



p' (x) = Iix^ - 6x2 + 1 
X = -1 



-6 



X = 2 

0 ' 



0 1 



10 -10 



-10 10 -9 , p.(-i) ^ -9 
i 

-6-0 I ' I 



p(2) -9 



Forquartic polynomials: 10 DIM A(5) 

20 MAT INPUt'^A 

.30 INPUT X 

^'kp LET P = A(5) 

50 FOR Ka'i TO 1 STEP -1 ' 

60 LET P = P*X +A(K) 

70 NEXT K 

80 PRINT X,P 

90 END 



16 




* p(t) 



p(-T) ^ 



10 110 1110 



11 111 nil p(io) = 1 , 



1 



10 , 20 



10 21 • (10101) = (21) 



s 1 



1^ 



1 3 7 15 t (nil) = (15) - 

• • 2. 10 

\ 



8 T6 



^ 8 16 (lOOOOL^ = (16) 

.*?.2- , 10 



3' 



56 ^80 38^0 30728 



60"^^.i♦80^ 38^1 *3073l (7^013) = 'boTii) 

• • • ' - 1 0 



5 \ 



' 3 



12 



'""Tl^ 



20 120 600 



^2 129 ^^^r ^<^^ ^ ' ' 



^3 129 m^^" ^(^^^^^^^^(^sa)- 



2i* 120 601 



3x^ + 2x + 6 



X - 2 I 3x^ - i*x2 + 2x - 3 



3x^ - 6x^ 



,2x^ + 2x 

2x^ - ^x 
- " ft — 

6x 
6x 



(^^01) , = (6-01^' ' . Xj-tf ? 



3x ^- i*x^ + 2x - 3 = (x -'2)(3x^ + 2x + 6) + 9 



10. 



2x^+ 2x2 - 



ix^ 
2^ 



X + 3 I x^ + x^ 



<^'^+ 2x3, 
x^ + x^ -ix^+lx^ 



X +2 



19* 



1 I* 1 3. 1 2 3 
-TX - yxHyX^ -l-x 

"1 
1 



3 1 2 . 



x" + x" 



X + 2x^ - X + 2 = (2x^ + 2x^_ - X + 3)<|x' 



x.l) 



ii 



= 1? 
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11. For n i it, m = 2: 10 DIM-A(5) .U(3) .^(3) .R;(2) 
^ 20 MAT INPUT A 

30 REM THE ROUTJt^E DESTROYS A 
^♦0 MAT INPUT U 
^ . 50 FOR K = 3 TO 1 STEP*-1 

60 LET (1(K) = A(K+2)/U(3) 
I 70 FOR K+! TO K STEP -1 

: 80' LET A(J) = A(J) - (|(K)AU(J-K+1) 

. • ' 90 NEXT J 
100 NEXT K 
^ - 110 LET ,R(1) = A(1) 

^ 120 LET R(2) = A(2) 

130 MAT PRINT Q 
1^0 MAT PRINt R 
150 END 



12. Suppose there were more than one pair. 

That is suppose ;)(x) = q (x)u(x)-4-r (x) 'where either r(x) = 0 

or deg r(x> < deg u(x) 
p(x) = q*(x)u(x) + rA(x) where either rA(x) = 0 
^ ' or deg r-^(x) < deg u(x). 

(q(x) - fli^(x.))u(x) - r(x) - r*{x) ^ 

.»» u(x) divides r(x) = r^(x) and thus deg u(x) < deg[r(x),- r*(x)] 

or (r(x) - r*U)] = 0. 



But 0 < deg [r(x) - rHx)] < deg u(x) so we must h^e r(x) 7 r^(x) 
or r(x) = rA(x)r 

*«. 

=*■ q(x)u(x):= q*(x)u(x) * 



And hence, q(x) = q"(x). 



Therefore there was real|y only one pair, i.e., th«t q (x) and r(x) 
are unique. t 



ERIC . * 



15. For n * k: 10 DIM A(5),P(5) 
\. 20 MAT INPUT A 



30 INPUT X . 
'fO FOR K-1 TO 5 ^ 
50 l,ET P(K) - A(5) 
60 NEXT K 

70 FOR K « 5 TO 1 STEP -1 

80 LET P(l) - P(l)*x + A(K) 

90 FOR J « 2 TO K 

100 LET P(J) » P(J)*X + P(J-I) 

110 NEXT J f 

120 NEXT K ^ , 

130 PRINT X 

1^*0 MAT PRINT P 

150 EHD 



n 1 



umap 



\ *t UNIT 264 J 



\ 



ALGORITHMS FOR FINDING ZEROS OF FUNCTIONS 



by. Werner C. Rheinboldt 




ERIC 



COMPUTER SCIENCE 



edC/umap /SSchapel st/^wton.mass 0S160 



AL^ORI'IHMS 'I'OF^ MXDING ZI.ROS OJ^ I^^^NCTIONS 
by 

Werner C. Rheinboldt 
Department of Mathematics and Statistics 
t University of Pittsburgh v 

Pittsburgh, PA 15260 



Table of Contents • 

1. SOME MODEL PROBLEMS *. . . 1 

2. EXISTENCE QUESTIONS ^ li 

3. THE BISECTION METHOD 8 

SOME LINEARISATION METHODS . • ^. tt 

5. RATES OF CONVERGENCE' . . 1i* 

6. A PRACTICAL ALGORITHM -.18 

7. REFERENCES . 

* . ',*->' 

8. ANSWERS TO EXERCISES ^ * . ./ . . . - 

* ^ ^ ; , * . ^ \ - , 



i 



Intermodular Deacription Sheet : oyAP Unit 264 

Title : ALGORITHMS FOR FINDING ZEROS OF FUNCTIONS 

Author : Werner C. Rheinboldt 

^Department of Mathematics and Statistics 
University of Pittsburg 
Pittsbyrg, FN 15260 

Review Stage/Date : (II , 12/15/78 

Classification: COMP SCI - . 

— I 

Suggested Support Materials : 

References: See Section 7 of te\t. * ^ 

' ~ • ' 

Prerequi s i te Sk1l I s : .j.^ 

1. Be familiar with the Mean Value Theorem. ' < . 

2. Be* famDlac w| the Intermediate Value Theorem. 

3. Be able to differentiate elementary functions. 

^. Be familiar with ojaking estimates usino absolute value notation 

Output Skil Is : - 

t. Understanding of standard bisect iph, secant, and Newton methods 
of root finding and appreciation of their 1 imi tat ions and 
st rong. pofnts . Introduction to more recent "root finding 
methods. • ^ 

Other tReigtgd'*1rnits ; 

Horaer's Scheme and Related Algorithms (Unit 263) 



/ 

i 



ERJC 



© 1978 EDC/Project UMAP 
All Rights Reserved. 



I 



• MODULES AND MONOGRAPHS IN UNDERGRADUATE 
MATHEMATICS AND ITS APPLICATIONS PROJECT (UMAP) 

The goal of UMAP is to develop, through a conwnunity of users 
and developers, a system of instructional modules in undergraduate 
mathematics and its applications which may be used to supplement 
existing courses and' from which compJete courses may eventually be 
bui It. * . . ' . 

The Project is guided by a Nat i^onaT Steering Corrnii ttee*of 
mathematicians, scientists and educators. UMAP is fUnded by a 
grant from the National Science Foundation to Education 
Development Center, Inc., a publicly supported, nonprofit 
corporation engaged in educational research in the U.S. and - 
abroad. ^ . 



PROJECT STAFF 

Ross L\ Finney 
Solomon Garfunkel 

Fel ic5.a DeMay Wei tzel 
Barbara KelczewskI 
Edwina R. Michener' 
Dianne. Lai ly 
. Paula M. SanXi 1 lo - 
Carol Forray 
Zachary ZevI tas 

NATIONAL STEERfNG COMMITTER 

W.T. Martin 
Steven J. -Brams 
Llayron Clarkson 
' Erhest J, Henley 
- Donald A. Larson 
Wi ri i am F. Lucas 
R> Duncan Luce 
George Mi I ler 
^Frederick Hostel ler 
Walter E. Sears 
^ George Springer 

Arnold A. Strassenburg 
Alfred B. Wil Icox 



Director 

Associate Director/Consortium 

Coordinator 
Associate' Di rector for Administration 
Coordinator for Materials Production 
Editorial Consultant 
Project Secretary 
Financial Assistant/Secretary 
Technical Typist/Production Assistant 
Order Processor ' 



MIT (Chariman) 

New York University 

Texas Southern University 

University of Houston 

SUNY at Buffalo 

Cornell University V^^.>^ 

Harvard University . 

Nassau Community College 

Harvard University 

University of Michigan Pre&s 

Indi ana Un ivers i ty 

SUNY at Stony Brook 

Mathematical' Association of America 



The' Project would like to thank Charles Votaw and Douglas F. 
Hale fer their reviews, and all others who assisted in the 
projiuction of this unit. " ~' 

This material was prepared with the support of National 
Science Foundation Grant. Np. SED76-19615 A02. -Recommendations 
ej<pres§ed are those of the author and^do not necessarily reflect' 
the views of the NSF, nor of the National Steering Committee. 



ALGORITHMS FOR FINDING ZEROS OF FUNCTIONS 

by ' 

Werner C. Rheinboldt 
Department of Mathematics and Statif^tics 
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1. SOME MODEL PROBLEMS 

^ Let f be some real function of a real variable 
X. want tjo find a real solution (zero, root) x 
of the equation 

(1) f(x) = 0.. . 

Only in a few cases, such as for linear or quadratic 
functions f, are there any explicit formulas for such 
a solution. Hence we will have to be 'satisfied with 
computing x approximately. ^ 

Nonlinear equations arise frequently in 
applications. For later pse we give here a few simple 
examples. „ * ' * ^ 

Van der Waal*s equation of state for an imperfect 
gas has the form 

(2) . (p + ^^)(v - b) = RT. 

v 

Here p [atm] is the. pressure, v Qiters/mole"] the molal 
voluiTie (volume/mass)^ T[Pk3 the absolute temperature, 
R F d-. 0820^54 [liter ~Utm/moTe^ .the "gas"conrta:ht"7 and 
a [liter.^ atm/mole^j , b [liter/mole] constants dependent 
on the particular gas. For instance, for carbon dioxide 
we have a « 3\S92, b = 0.04267 and for helium 



a ^ 0.03412, b = 0.02370. 

For given values of p, T, a, b we want to" compute 

the corresponding value(s) of v for which Equation (2) 

holds. This is problem of the form of Equation (11. 

2 ' 

More specifically, after multiplying by v , tha.^desired 
values a-re the solut ion^^'^^^f^Jfie^ubic polynomial in v 

(3) pv^ - (pb .+ RT)v^ + av - ab = 0. 

As another example consider the motion of a 

particle of mass m which is attracted to a fixed center 

2 

0 by a Newtonian force pm/r with constant u > 0. 
Kepler's first law then states that the particle moves 
oir a conic section with eccentricity e with one focus 
at 0. . Thus for 0 < e < 1 the orbit is an ellipse, for 
' e = 1 a p4.rabola, and for e > 1 a branch of a hyperbola. 




I 

Figure 1, 

y 



More specif i<j[ally, let P be the pericenter, that is the 
point on the orbit closest to 0, and introduce the polar 
co6rdinates (r,<>) with center at 0 and the direction of 
OP 'as the x-axis . Then for 'e > 0, e ^ 1 the orbit is 
given by JO , * . 



2 



Now let T be the^time when the particle is at the 
pericenter, then it| ^position at time t is determined 

by the ICepler equations 

-* 

(a) n(t - T) = u - e sinu, if 0 < e < 1 H 
} (b)' n(t - T) =- esinhu • u, if e > 1. 

The va?*iable u is called the eccentric anomaly; it 
relates to r by the equations 

* r = a(l - e cos u) , if 0 < e < 1 * 
^ ^ r = a(ecoshu - 1), if e > 1. 

The parameter n is the mean motion, that is, la the 
case of an ellipse n = Z-n/p where p is the period. For 
given^a, e, n, T the problem of determining the position 
of the particle at time t now requires the solution u of 
the corresponding Equation (5). Then r can be found from 
Equations (6) and <J> from Equation (4). 

Some values of the parameters for the case of the 
earth's orbit around the sun are a = 1.000, e = .017, 
n = .017?0, and T = Jan* 1, 1900. In th^s elliptic case 
Equation (5) (a) is unchanged if we add or subtract a 
multiple of Zir from i .= n(t - T) and u. Hence we may 
^ always reduce the left side such that -n < ^ ± n. 

Exe rcl ses " . ' 

1* Let the 'function g be continuous on the closed interval 
""^a X £ b and di fferentiable on a < x x b.. Then the Mean 

Value Theorem ensures the existence of at least one value ' ' 
x'^ such that* ^ , * « 

Thus to find x* we need to solve some nonlinear equation of 
the form of Equation (l). ^Write down this equation in the 
• following cases: 

' ' 

. (a) ,g(x) - (x + J)e"^, a » -1, b « 0. 
^ (b) g(x) « x^^, a ^« 0, b « 20^^^^^i 1. 170a. 

FRIG oi. 
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2. For a continuous function g on the interval 3 x <^ b there 
« exists at least one value x such that 



g(x)dx = (b - a) g(x ) ,*a» < x < b. 



This is the integral mean value theorem. Write down the 
resulting equation for yr^ in the cases 



(a) 



dx 



X log X ' 



(b) 



I 

' (x + e^) dx. 



2. . EXISTENCE QUESTIONS 

-'Before we look at methods 'for solving a given 
Equation ^1) , "^it is important to realize that there' may 
be no solution at all or there may e^ist any number of 
them. The examples in Table 1 illustrate some of the 
possibilities . 

TABLE 1 . • 



•f(x) 



No., of zeros 



Zeros 



rx - 1 



x^-i 



tan X 

,2 ■ 



(xJ-1) ^ for. X ^ -1 
0 ig^'l li X < 1 



(x-1)^ for X 2. 1 



one 
two 

countably many 
a continuum 



X = 2 



1, X. 



-1 



Xj^ = kn, k=0, +1, +2, 



any x with -1 1 x* <^ +1 



This indicated that it is advisable to begin an 
investigation of a particular equation by 'localizing a q • 



suitable interval a < x < b m which the desired solution 
exists. . A simple approach for the construction of such an 
interval is to plot the graph of f(x). For instance in the 
case of the simple Kepler equation 



(7) 



f(x) 



-x - 1 



0 



this may provide the results shown in Table 2 and -Figure 2. 

TABLE 2' r 



X 


f(r) 




0 ' 




-1 


0000 




0.5 




-0 


'7397 




1.0 




-0 


4207 




1.5 




0 


1252.10"^ 




2.0 




0 


54 54-X^ 




. 2.5 




1 


.2008 




3.0 




1 


.9294 




figure 2. 



The drawing indicates that a z,eio of f should be contained 
in the interval 1 <_ x <^ 1.5. 

The theoretical foundation for this conclusion is the 
following basic theorem: 

Intermediat-e Value Theorem:. Let f be a continuous ^ 
^ 

function on the interval a <^ x <_ b , a < b , If 

(8) ~ sign T(a) ^ sign f(b) ^ - 

then fix) = 0 has a solution in the interval 
a ^ X , _< b . _ ^ 

The proof of this fheorem is not entirely simple 

and uses some fundamental properties of the real number 

system. But intuitively the result is clear. A continuous 

function might be charac reri zed as a function wi,th a graph 

that can be drawn without lifting the pen fronj the paper. ^ 

The condition (8> implies that at two endpoin«^' of the 

interval the function values are on opposite sides of the 

X-axis; hence when drawing the graph of f we need to cross 

* ' 'ft 

the X-axis somewhere in that interval. * 




Figure 3. O^^ 





about tlie "nwber o-f so\\^ii>ri^^^^ ' : ^^.'^ >: r; ' ^ " -i 

a <^ X '<_ "bv fact, tlio.ve 'jnayy-ir^ 

Vo^t s . ^ For Jn4 tiiv^.'e'f.. in .the- i^it'i^fWi . -.^ -'S4% 1 V' " 

(9) - ' f{x> ^' x5 .\ ^ -jv/-/ ^ 

has for a ^= - l/2:.,a /:^*Xagie .root /i-r/^.-^^^^ two ^.y-': - 

roots' X = O and X s 'lv^-Z^'-S^T^^^^f i'lra^^^^^ A/34, -r^^ 

exactly three roots. ' ^ - " ^ 

: - ' - ' r ,^ - * 
The c>ase a = 0 .is exceptioTvaj-; Gejie rally,, a^robt x 

of tfie Equation (1.) is. called a root oiT; multipli'ci^ty m/ff 

for X near x, the function f cSfn be written in the 

(10) " ^f(x).= (X - x*)^^(x)-^ 

with some, function g that is continuous near x- and 
satisfies g(x*) i- Q. In our case, Equation (9) has for 
a = 0 the ¥orm f(x) = x^(x •'Jj) which shows that x* = *0 
IS a root of multiplicity two while x «= has multiplicity 
one^. Thus counting multiplicities we realiy have three 
roots. Generally, the followin^g result holds: 

Theorem : Under the conditions of the intermed'iate 
value theorem, the interval a <^ x < . b contains either 
infinitely many solutions or finitely many solytions 
for which the^%um of their multiplicities is an odd 
number ' - ^ — . ^ - - - ^ 



Or- 




" 1 ~~ tfefermi-fie the solutions of .the f un.pt loji' 



fpn^x = 0 



_X .5 i^h' — fdf^^ < X < 1 



rn- .thei^. Interval 0, <,^,.< i.'-'' 



Fo;- any 0 e ' < J and't^ -<'2, < tt det^rfnlne- an interval which 
contains a root of the .Kepler , equation £ » u - e sin i/. 

-3;^' ^howTt hat Jqr t K^' So' *c a lied c i* i t i c a 1 va 1 ue s 

a , 8a ' ' 



- c 



J - "^c = 



27^:^;-^ " 27Rb 

van der Waal* .Equat ion (3') haS only one root v^ = 3b of multiplicity 
thV^,.' Show that the constants a, b, R can be expressed in terms of 
the cr i t Leaf" va lues as 

^ , 2 ^ 1 „ 8 ^c^'c 
a = 3p,v^ , b = ^v^ , R = y ^ . 



Use this to show that with the d rmens r on- 1 ess variables p = p/p , 
» * c 

V = v/v , T = T/T the -ecfuat ion assumes the form ' 
- c ' c - 



(p + 4") (3O 

V J 



r) = 8T. 



3. THE BISECTION METHOD 



Suppose that an interval a <^ x < b has been 

found >vhfere the conditions of the intermediate value , 

fheor-em are satisfied. Then we know that there is at 
t 

least one solution x of Equation (1) between a *and b. 
, Fpr'the midpoint m = a + (b" - a)/2 we test now the 

condition sign f(m) f sign'f(a). If* -it holds then the 

intermediate vaKie theorem guarantees that there is- a 
^ root .^n the interval, a < x < m, otherwise , we have '* 

s*ig;n fCm) * sign f(a^ f sign f(b) and hence there must 



Of; 



be a root in m < x < b. In either case, tixe 1 ength 
of the interval has been halved. ^'B^^epeat ing 
^prQcad^we can decrease the "interval - length below a 
prescribed tolerance and henge approximate a root of f 
arbitrarily clQsely. 

In an informal programming language this algorithm 
may be written rn the following form. 

1. Input (a, b, kmax, tol}; 

2. k^:= 0;! . ^ 

3. If (sign f(a) = sign f(b)) then error return 1: 

"Wrong interval 

Print {k, ^, b}; . 
5- U Jt>-a| <^ tol then normal return; 

6. "k := k+1; 

7. m := a + (b-a)/2; 

8- 11 (sign f(a) ^ sign f(m)) then b=m else a=m; 

9. U k < kmax then go to 4 else err^r return 2 
- , Vkmax exceeded**; 

Step 2 has been included to verify that at all times the 
basic condition (8) is satisfied. If it holds for the 
- input interval then theoretically it will remain valid 
for all subsequent intervals. But in practice this may 
well not be the case due to round-off errors. All 
iterative methods should include a count k of the number 
of steps taken and use , it to stop the process when a 
given maximum count kmax has been exceeded. This is ^ 
done here in step 9. 

Table 3 shows the results of the first five steps 

when the algorithm is applied to the Kepler Equation 

(J) on the interval 1 < x < 2. 

— — «i 

4|P^ 



Table 3 



k 


a - 


b 


b-a 




f(b) 


•0 


1 


2 


1 






n 


d c; 1 

J J ^4 




1 


1 


1 . 5 


0 . 5 


-0 


.4 207 


% 


IZSJ ' 10' 


2 


2 


1.25 


1.5 


0.25 


-0 


.'2245 


0, 


1252-10' 


2 


3 


1. 375 


1.5 


0.125 


-0 


.1154 


0, 


1252-10' 


2 


4 


1.4375 


1.5 


0.062 5 


-0 


.05806 


0 


1252-10' 


2 


5 


1.46875 


1.5 


0.08125 


-0 


.02865 


0 


1 252^-10' 


2 



We shall see later that, with eight digit accuracy, the 
exact root is x = 1.4987011. Obviously, otfr algorithm 
is not particularly fast. 

• The interval decreases at each step by a factor of 
two. Hence the kth interval has the length (b - a)/2 . 
If the tolerance is, say, 10 then we require that 



or 

that is 
(11) \' 



(b-a)2 <i 10 



2'^ > (b-a) 10"-, 
k,> log2r(b • a)10^J. 



lO the example of Table 3 this means that we need 
k = 23 to obtain seven digi t 'accuracy . 



Exercises 
1. 
2. 



Draw a flow chart for the bisection algorithm. 

If a programmable calculator or computer is available, 
implement the algorithm for ^t,he general Kepler Equations 
(5) with given = n(t ~ T) <^ tt and e > 0. 

3. For the cubic polynomial f(x) = x - 2 x + 2 determine an 
interval containing (only) root and apply the bisection 
algorithm to approximate the root to four digit accuracy. 
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4. SOME LINEARIZATION METHODS 



Jn order to overcome the relatively slow convergence 
of the 'brsection method, we turn now to a deferent 
principle for computing solutioi^s of Equation (1). It is 
based on the idea of replacing the function f(x)^by a 
succession of linear functions g].(x) = aj^x + h^, k^,2,..., 
such that their zeros . -b^^/aj^ , k = 1,2,.... approximate the 
desired solution x of Equation (1). 

A linear function is determine^l by Its^values at two 
distinct points- .Suppose that we'ar'bat the kth step of 
our process ^nd that the aj^proximat ions x^, x^, of 
X* have been computed already, 1 f, k > \ and Xj^ f ^Y-l 
then we may construct the l,inear functior] g}.(x) which 
agrees with f(x) atf Xj^ and Xj^^ ^ namely 

fcx.) - ax, .) 

.(12) '&kW = f^y- X, ^ X,,, ' ' ^k)' 



For f(Xj^) °^ ^^^^ function, (that is 

i • • ^k ' ^k-1 ^ • . • 

(13) x^V^ = x^ - f(x^) - f(x^^.^) ^^H^' 

is taken* ^s the next approximation of x and in this way > 
. the process is continued. This is called" the secant method' . 
' Clearly the process will fal^l when two equal function values 
•ana encountered: But if this does not happen then we 
shall see that the nfethod is considerably faster than the 
bisection method, . * ► 

A. linear, function is also defined by its function^ 
^ value and i^s slope at a given point. Hence, suppose that , 

at Xj^, k > 0, we are able- to compufe no,t only f(Xj^) but 
- also the value. of th^def i vtfti ve f*(x,)., * Then we can replac 
tne secant line. (12) by the tangent riae . 

(14) g^(x) = f(x^) + f'(^k) (x-Xj^)." 
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Jf f'CXj.) 0 then the zero of Equation (14) is 
(IS) . x,,j = X, ^— 



which becomes our next approximation of x . The resulting 
process is called Newton's method . It will fail whenever 
a zero derivative value is encountered, but otherwise it 
turns out to be even faster than the secant method. 

As an example, we consider the computation of the 

positive square root of some positive number a > 0. In 

other words, we wish to find the positive root of f(x) = 
2 

X - a. In that case, we have 



f(x^) - f(x^.j) - xj.j- 



and hence the secant method assumes the, form 
On the other hand, Newton's method becomes 

^''^ '^k.i = -k -'^ = 

For a = 10 we start (16) .with x^ = 11, x^ = 10 and 
(17) with = 10. The resulting first few steps are 
^given in Table 4 . ' '* 
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Table 4 



Secant Method 



Newton 's Method 



^0 
1 
2 
5 
4 
5 
6 
7- 
8 



11 
10 

5.71429 
4.27273 
3. '44^,0 3 
3.'20510 
3.164^0 2 
3.162^9 
3M6228. 



Ill 

90 

-.22.6531. 
8.25620 
i'. 87515 
0.2S9849' 
0.110'2U' 10" 
0.7064-10"'^ 
/ 0.2-10'^ 



10 ' 
5.S «' 
3.65909 
3. 19601 
I 3.16246 
'3. 16228 



.90 ' 
20 . 25 
3. 38895 
0.214448 
,112557'10" 
0.2-10"'' , 



Exercises. 
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K Apply the secant method to the Kepler 'Equat ion (7) 
starting with Xq » 2, x^ « 1.5. 

2. Use Newton's roetho^ starting from Xq»-1 to solve the 
equation specified by problem 1(a) of Section I. ^ 

*3. If a prograninable calculator or computer is available 

Implement Newton's method for the computation of the" ^ 
square root of any positive number a. Use Xq a as 
starting point. 

'Tor polynomial etjuations the* value of the function and 
its derivative at any given point may be computed 
^^>lmultaneously by cleans of Horner's Scheme.* Draw a flow 
"chart of the resulting process. If a "programmable 

calculator or computer Is available implement the method 
• for cubic ar^d quartic polynomials and test it orl severaf 
' ^equations, sUch as the polynomial (3). 

or explanation of Horner's Scheme, see UMAP #263. ' 



< 5. RATES OF CONVI-RGHNCH 

The bisection method generates a seq^ience of intervals 
_< X bj., k = 0, 1, . which contain the desired 

rooi X*. Any point in the kth 'interval ipay be considered 
as the kth approximation of x ; for the moment ^let us 
consider the midpoint m^^ = a^^ + (bj, - for tha4t 

purpose^ Since at each step the interval is halved we 
/then have the obvious relatibn 

(18) hj^^j -^i <_ \\m^ -x*|, k = 0, 1, .... 

In other words, the errors converge to zero at least as 

* • * k 

fast as the geometric sequence -x |/2 , k = 0, 1, .... 

I Now s.uppose that Newton's method is used and produces 
a sequence of points Xq , x^,. X2, ... which converge to the 
solution X of .g^g^tion (1) . >loreover, assume that " 
the Xj^^are all contained in somp interval a <-x < b 

where 

(19) (i) |£'(x)| >*a > 0', /(ii)- |f'Ux)| < 6, 

for a X <^ b. ^ 

Then it follows^y Taylor's formula that with certain 
Cj^ in our interval 

0 = f(x*) ='f(xj^) + f(xj^)(x* - Xj^) ^ ^f"(C^)(/ - Xj^)2 

By tHe definition (15) of Newton's method the term in the 
■s^ua're bracket is ze,ro whence 

0,= f«(Xj^)(x*-- Xj^^p + ■?f"('^k^^''* " ""k^^ 

or, because of condition (^9), • 

(20) |x - X. I = f- ^ Ix -.x^l^- < II |x, - X |2 

k+1 I |f(xp| ^ Z a 



In other words, in this case, the error is squared at 
each step (up to some factor). Thus if, say, the k 
error is proportional to lO"^^ then the »(k + l)st error 
is proportional to 10 . fhis is clearly seen in the 
third column of Table 4 where' For k ^ 3 the number of 
correct^ digits doubles with each step. 

For the secajit m'ethod a similar estimate can be 
derived.^ It turns out that if again Condition (19) 
holds and all x, remain'in a _< ^ b then* we have 
with same constant y > 0 ^ S 



(21) |x - X, J 1 Y '|x - Xj^r , k = 0, 1, . ..9 



* I J * 1 1 

t = |(1 + /T ) ='1.6180. 

1^ Thus the errors tend to zero somewhat more slowly than 
in the case of Newton's metKpd but certainly- faster than 
with the bisection method.^ | ' 

These results are somewhat deceptive. First of^ll., 
the condition (19) (i) is fairly strong since it implies 
that there ,is only one root in -the interval a <^ x <^ b 
and this root must have multiplic^y one. In fact, by 
the| mean-value theorem we have ^ • 

|f(x)| ^= |f(x) - f(x*)| = |f»(0(x-x*)| > a 1x^- x*| , 

* a < X • < b«r • 

Hence f(x) = 0 for any a _< x ^ b implies that x = x . 
Moreover, if f can be written in the form of Equation 
(10) with g(x ) ^ 0 then the left ^ide of 

|x ' x*r-^|g(x)| > a > 0 \ 

tends to zero, as x goes to x unless ni = 1, 

For zeros of multiplicity gre,ater than one, Newlr^n's 
method^ indeed} converges much more slowly.. For exam^Le, 
in the simple case 

' {22y f(x) = (X -^1)"^ = 0, X* = 1 

^ "jw.ton^s method has. the' form 15 



(^k -^^"^ m . 1 

/ , s m- i 
m(x, -1) 



whence 

a 



'k+1 " m '^"k 



^x, 1). 



In other words, foi: .m = 2 the convergence is here as 
slow as that of the bisection method, and for m > 2 it 
is even slower. The secant method shows a similar 
behavior; in fact, not only the rate of convergence 
deteriorates but arbitrarily close to the root we may 
encounter f(Xj^) = ^^^k-l^ which case the method 
fai'/s completely . ' , ' 

There are further problems with the function (22) . 
In fact, we see that |f(x)| _< e implies tliat 

1 - ei^/"^ 1 X 1 1 ^ e^/"^. 

-iience, say; for m = 10 we have |f(x)|'"*^10 for 
.75 < X <i 1.25. In other words, with six digit 
accuracy any point in this interval can be called a 
zero of f. Unless higher accuracy is used any iterative- 
process entering this uncertainty interval will, by 
necessity, show erratic behavior. A root of this kind 
is called ill-conditioned. It turns out that also roots 
of multiplicity one may be ill-conditioned.' * 

Even if the conditions (19) hold the estimate (20) 
for Newton *s method may be very misleading. Consider for 
example the equation 

(23) x^^ -'-1 = 0? - ' 

we encountered in Exercise 1(b) of Section 1. Herff^ 
clearly f'(x) > 0 for x > 6 and for any interval ^ 
0 < a <_ X b containing x ^ 1 the estimate (20) 
holds. But the factor 6/2a will be very large unless 
a "and b are very close to one. This reflects difficulties 
w^ith Newton's method, and in fact for Xq = 1/2 we have 
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= 13,797.53 which is certainly a much worse approximation 
of X* than Xq. The subsequent iterates decrease 
monotonically, but very slowly, to one. Only very close 
to one the expected rapid conver^ice sets in. 



Exerci ses 



and 



1. Show that when Newton's method is used for solving 
- a » a > 0, starting from any Xq > 0, 
Xq /a , the iteraites satisfy 

^ > X2 > ... > > x^^, >/a , kl 1 

and 

2., Apply the secant method to th6 equation 

f(x) « x^.- a « 6, a > 0, starting from Xq > x^ > /a". 
Show that ^ 

- (x,-/iy(v, -vi) 



3. Apply Newton*s method to Equation (23) 'starting from 
Xq « 1/2. Show that 



18 1 : Jl ^ 

. How many Iteration steps are needed to reach |x. - l| ± 1.17 



6. A PRACTICAL ALGORITHM , 

The results of the previous sections show that none 
of the methods discussed here is entirely satisfactory. 
The bisection method is 'fairly reliable but slow, the 
Newton and secant method are both much more rapid in 
certain cases but show unreliable behavior in many others. 

We discuss now an algorithm which combines the 
bisection and secant .nethods to bring out their best 
^atures. ^ Again we work with a sequence of intervals 
0^ decreasing length for which the intermediate value 
theorem holds. If, say, a _< x <_ b is the kth interval, 
then we set 

Xj^ = a, yj^ - b if |f(a),l- ^ |f(bO| 

Thus Xy, may be considered the current best approximation 
of the root in the kth interval., 

A step of the algorithm now consists in determining 
a new point between Xy, and y^^, called w f<j^ the' moment. 



which will become either x 
•introduce the point 



k+1 



or y 



k+1* 



For this we 




x^^ ^ if k > 1 and y^^ 



otherwise 



any consider first the secant step 
^2^) • ^ = ^k - f(xv) ■ f(zj 



provided it gives a better result than the bisection steo 



(yj, - x,)/2. 



In other words, since Xj^ is the current best ajyproxjmation , 
s has to be between x^ and m/ -At the same tiipe, since X|^ 



is not yet within a given tolerance of the root, s should 
differ from Xj^ at least by that tolerance. 

Before we discuss the choice of s or m, a few words 
about the definition of Zy, may be useful. The normally ^ 
eicpected choice would be = yj^. Here^ represents a 
secant step based on the two current endpoints of the 
interval where the function values have^ different signs. 
This is called a regula-falsi step. Unless 'rourtd-off 
interferes," such steps do not lead out- of the interval 
and, since they have'no subtractive -cancellation problem 
In the denominator of (24), they are generally rather 
stable. But in situations, such as that shown in Figure 
4, regular falsi steps may .give very poor improve^nts 
, of the interval. For this reason, we use in the case qf 
^ ^k-1' ^ - ^ secant step based on ^^^^ and 
Such a step may lead' entirely out of the interval and 
hence has' to be carefully controlled, but it certainly^ 
guarantees that there will be ivo long sequence of small 
.steps of the type shown in the figure. 



V 




Figure A. 



In order to test convergence we use the tolerance 
.funetion 

i 

tol(x) = £|-x| +6 
where g > 0, 6 >^ 0 are given constants with 



e 



^ a > 0. ' For-G =' 0 the condition |x - x | < tol^ (x) 



requires the absolute error |x x | to be below 6 
* while, for a = 0 it forces the relative error 
gl^^-x - x*j/|x| to be less than g. ' ^ 



it 



With this we set now w = s if s is between 
t = Xj^ + sign ( yj^ - x^) tol ( Xj^) and m, and w = t if s 
is between Xj^ and t. In all other-cases, w = m is 
chosen, that is, we take a bisection step. 

Thus in either case we have settled on a value of 
w. If sign f(-w) / sign Hv^) then the interval between 
w and Yy^ is our new mterv.al, otherwise the interval 
between Xj^ and w i's chOsen. -This completes one step -of 
the process . 

We terminate the algorithm if the length of the 
interval between Xj^ and m is less than tol ( Xj^) , that^.is, 
if \y^ - xj^t 1 2 tol (\)% This fits with-our choice of 
the minimum step, w = t wll«p s is between and t. In 
fact* if we have not yet Converged then 
^ |m - Xj.1 > toKkj^) = I't - Xj^^ . ' • . 

and thus^also in tHi§ case w is always between x^^ and m. 

•/pr the implementation of this process we have to 
take care;,that the division in^fthe secant step (24) does 
n^t produce overflow or jtnderf low. For this we compute 
the numerator and denominator , 

separately and then test ^ * 



j\y^ x^lkl V [p|-l |q| tol(x^) 
to determine whetther s will be between t and m. 

The overall algprithm can noW be formulated as 



follows . 



1. Input {x,y,e,6,kmax} ; * . ' 

z<^; ^ ^ ^ * ' 

3. *k:=0;' * ' ' , 

4' » If (|f (x) |>|lP(y) I) then z:^x-, x:=y; y: = 2; 

5. » Print ^rk,'x,y} ; ' ^ 

6. ; tol:=eUl +6;' ^ .^'3, ' 



7. U (|x-y| < 2 tol) then normal return; 

if ^ kmax) then error return "k max exceeded" 
9. k:«k+l • 

10. ,p;=(x-2)f(xj; 

11. q:=fu) - f(x); 

12. U ip < 0) then p:=-p; q:=-q; 

13. 2 : =x; 

14. if (p < jq(tol) then x:=x + signt>'-x)tol 

^ » else ij; (p < b(y-x)q) then x:=x + p/q 

' . ' • • else x : ='x + ( )' - x ) / 2 ; 

' 15. (sign f(x) - sign f(y)) then y:=2;- 

16. Go to 4^; 

, As an example we gave in Table S the solutior> of 
Equation (23) with € = 0, <S = 0.5 x 10'^, and. the starting 
interval x'= 0.5, y - 2.0. the straight bisection method 



« uses 


in this case 


about 21 steps. 














Table 5 










-A 


Interval • 




F(x) 








0.500000 


2.00000 


-0. 


^99998 








o^soooop 


1.250000 


-0. 


999998 








0.^10811 


1,250000 


-0. 


999997 








0.880406 


1.250000 


-0. 


911084 








0.880406 


. 1.06520 


-0. 


911084 






> 


0.932504 


1.06520 


* -0. 


734927 








0.998854 


1.06520 


-0. 


215500 


•10 


-1 




0.998854 
'0.999991 


3^.00086 
1.00086 


-0. 
-0. 


215500 
168309 


•10 
•10 


-1 

-3 




1.00000 - * 


1 .0008^^^^ 


-0. 


127405 


•10 


-S 




1.0.0000 


1.00000 




127405 


•10 





The process docs not always perform that well. < 
In fact, for instance, for. Equation (22)>'ith m = 19 
and X = 0 , y = 10 , e = 0 , 6 = 0. S • 10"^^ it requires 
several thousand steps. (Fven for 6 = O.S-IO'*** a total 
of 121 stcp^ are taken). The reason is that for a long 
tunc the algorithm uses only the minimal step^ of 
length tol ( Xj,) . ' * • 

Various remedies have been proposed for this 
problem. The easiest approach is to force periodically 
some bisection steps. For instance, we may simply add 
a test between steps 12. and 'l3. which forces a bisection 
step if k is a multiple of some fixed period M and then 
bypa3^ses s^ep 13.'' More sophisticated is a test every 
M s^eps which determines whether thfe interval at the 
beginning of the period has been reduced at least by 
the factor 2 corresponding to M bisectjon steps. We^ 
leave the details as an exercise. 



Exerc i se " . 

1. If a comput^ is Available implement the process describe^ 
in this sectlpn. Test it with the exanjple of Table 5, then 
use It to solve the equations of Section 1. • 

2.. Introduce hn your program a forced bisection s«tep when the 
i teration index k is a. mul t iple of some integer M > 1 . 
\ Apply the resulting process to Equ^ion (22) and experiment 
with different values of M. \/ 

3. Instead of the approach of Exefcise 2 introduce a^'ffgl^iod ic 
check comparing the actual reduction of the interval with 
the expected reduction by means of the bisection method. 
Compare the performance with that of the process developed 
f( in Exercise 2 above. 
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8. ANSWERS TO EXERCISES 



Section 1 



1 



(a) 
. (b) 

2^ (a) 
(b) 

Section 2 
1. X 



1 + xe 

. 1 



0. 



= log l2iL4 - 10.63289. 

X log X ^ log 2 



X + e 



1 

I 



1, 2r 



2.21828. 



, and X* 



0. 



2. For f (u) - I ' u f Q sin u we have f (0) = I and 

f(-TT) Ji + TT > 0, f(TT) = A TT < 0. ' 

Thus for 0 < Ji>- <• TT we may use the interval 0. £ u £ tt 
and for -it < A < P the interval -"t 1 u £ 0. . 

3. For*the critical calues p^ and T^ tlie van der Waal 
polynomial (3) reduces to 



0 



9bv^ + 27b^v - 27b^ 



(v-3b)' 



and thus has only one triple root v^ = 3b. 

From v = 3b it follows that b = v^/3 and hence from 
^ 2 2 

the formula for p^, we obtain a = 27p^(v^/3) = ^Pc^c* 

2 

Now the expression for T^ gives R = 24p^v^/ (9T^v^) = 
(8/3)p^v^/T^. ' By subs4:ituting these quantities int^o 
(2) we find that 



'P ,+ 



2, 



iPc ^c T 



which 'leads to the stated dimensionless equation after 
• multiplication by 3/(p^v^>. 

Section 3 ^ • ' 

3. '^f(-2) = -2, f(-l) = 3, X* i -1.7693 . 
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Section 4 ' 

1. The iterates are in sequence 

2.0,. 1.5, 1 .4988490, 1 .4987012, 1 .4987011 

and tjie last number is correct to eight digit^. 

2. The iterates are in sequence 

'1., -.68393972, -.57745448, -.56722974, -.56714330 
and the last number is correct to eight digits. 

3. A simple infornral program for this might look as . - 
follows : 

1. .Input {x, n, a^, a^, a^^^, kmax, big, tol}; 

2. k:=« 0; 

*3. Print {k,x}; 

4. k:= k+1; ' * . * 

5. p:= a^^; ^ # 

6. pprime:= p; 

Z* For i = n-1, n-2, 1 do 7.1 p:= p x x + a^^; 

7.2 pprime:= pprime ^ x + \ 

8. p:= p X X + a^r 

9. H (IpI >*'tpprime| x big) then error return 1: 

^ '^Excessive step" ; 

10. H IpI, < |pprime| ^ tol) then normal return; 

11. x:f X - p/prime; 

12. 1£ k < kmax then go to 3 > 

else error return 2: "kmax exceeded"; 

r 

Section 5 

1. By (17) we have 

and thus 
as well as_^ 
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Hence for all k ^ 0 we have the implications 

Xy. > 0, Xy, ^ /a --^ ^k+1 ^ 

^k > ^k ^ ^k*l 

from which the stated inequalities follow directly 
by induction. 



2. By (16) we have 

X 



Vk-i " ^ 



k.l x,^ - x,^.^ 



and thus 



^k.l - = x,^ ^\y^,^ tVk-1 - -^(^k " ^k-l^ 



3qp7^5^(x^- /^)(x^.i - /a) 



as well ajs 



^k " ^k+1 + X, , (^k " 



2 

^k " ^k-l'' 

Hence for k >^ 1 we have here the implication 

^k-l > /a. x,^ > =^ x,^ > x,^,^ > 
and the stated result follows directly by induction. 

3. Newton's method here has the form 

19 , io 19 - . 

- 1 ^ 18 1 

-k.i = -k ; = i9xj;^ '^^"^^T^ 

Even fo;- Xj^ = 1.1 the second term on the right 
is only of the order of 1/100, and hence until then 
the principal reduction comes principally from the 
first term. For x^ = 13,797.53 we have 

for k = 175. 

« 

Section 6 

3. An informal program incorporating such a periodic 
check might look as follows: 
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1. ^Inpu't {x, y, 6, kmax}; 

2. m:=0; 

3. length:=ty-xL; 

4. z:=y; 

5. k: = 0; * 

6. If (|f(x)t > {-f(y)|) then z:=x; x:=y, y:=r; /J ' 

7. Print {k, x, y}; ^ 

8. / tol:=ejx| + 6; ' ' ' j ^ . ' 

9. I_f (tx-y| < 2 tol*) then normal return; , ^ f 
10. 1^ (k *> kmax) then error re/turrt "kmax exceeded"; 



11. k 

12. p 

13. q 



= k+l 



:=(2-x)f(x) ^ , 

=f(z) - f(x) 

14. If (P < o) then p:=:p; q:='-q; 

15. z:=x; 

X6,-« mr=m+l ■? ' 

17. If (m > 4) then if (16|y-x| > lengtfi) 

then x:=x+Ji(y-x) ; go td^20;* 
else m:=0; length: = |y-x( 

if. (P 1 Iqltol) then x:=x + sign(y-x)tol 
else if. (p < ^(y-x)q) 

-4" ' then x:=x + p/q 

else .X : =x + %(y-x) 

19. If (sign f(x) = sign f(y)) then y:=z; 

20. Go to 6 . 



I 



^ -STUDENT FORM 1 
Request for Help 



Return toj 

EDG/UMAP 

55 Chapel St, 

Newton, MA 02160 



Student : If you have, trouble witiv*9 specific part of this unit, please fill, 
out this form and take it to your ihstructor f or^assistance. The information 
you give will help the author to tevise the unit-'. 

Your Name 




Ok 



ERIC- 



Section 



Pdf agraph 



OR 



Description of JDif f iculty : (Please be specific) 



Unit No. 



Model Exam 
Problem No. 

Text 
Problem No, 



'3 



Instructor ; Please indicate your resolution of the difficulty in this box. 
Corrected errors vin^materials* List corrections here: 



p\) ^^^^ student better explanation, example, or procedure than in unit, 
Give brief out^line of your addition here^• - 



3 



o 



Assisted student in acquiring general. learning and problem-solving 
skills (not using examples, from this Unit,) 

■I 



Instructor's Signature 



please use reverse if necessary*/ 



I z**^^^ Return to*: 

STUDENT FORM 2 .-EDC/UMAP 

• Unit Questionnaire - ^^^^^wa^5^, 

Newton, MA 02160 

Name ^ ' Unit No. .Date 



Institution {' C ourse No. 

Check the ch^iie for each question that comes closest to your personal opinion. 

1. How useful was the amount of detail in the unit ? 

Not. enough detail to understand the unit 

^Unit would have been clearer with more detail 




Appropriate amount of ddtail. 

Unit was occasionally top detailed, but t|f®^was not distracting 

Too much detail; I w^Skftkten distracted 

2. How helpful were the problem answers ? ' , 

Sample solutions were too brief; I could not do the intermediate steps 

^Sufficient, information was given to solve the problems 



_Sample solutions were ^oo detailed;.,I didn^t need them 



3. Except for fulfillilig the prerequisites, how much did you use other sources (for 
example, i-nstructor, "friends, or other books) in order to understand the unit? 

' A L^t , _^ ^Somewhat y A Little Not at ^11 

4. How long was this unit in comparison to the amount of time you generally spend on 
a lesson (lectj|re and homework assignment) in a typical math or science course? 

Much > Somewhat About Somewhat Much 
\^ ^Longer. ^Longer the Same ^Short^^^ ^Shorter 

5. Were any of the following parts of the unit confusing or distracting ? (Check 
as many as apply.) 

Prerequisites 

^Statement of skills and- concepts (objectives) 



^Paragraph headings 
Examples 

^Special Assistance Supplement (If present) 
Other, please explain 



6. Were any of the following parts'" of the unit particularly helpful? (Check as many 
as apply.) 
' P rerequisites'" ? 

Statement of skills and concepts (objectives) ' • ' 
^Examples 



J^roblems - , ' 

Paragra ph Ji/ adings , j\"Vf 

JTable 8T Contents ! ijf 

Special Assistance Supplement (If present) 

Other, please explain : 



Please describe anything in the unit that you did not particularly like. 



•Please describe anything that you found particularly helpful. (Please use the back'^of 
this sheet If you peed more spaceO 

ERIC . 



o 



■ Corrections 

UMAP Module 264 

20, line-12 ' ' 

> = (x,^ - z)f(Xj^), q = f(z) - f(x^) 

21 , lines 1 to 11 

•.V, 

7. * If (|x-y| 2 tol) then riiormal return;/" 

8. If (k >'kmax) then er^ror return "kmax exceeded"; 
j\9. k:=k+l 

10. p:=(x-z)f{x); 

11. ' q:=f(z) - f(x); . ^ 

12. If (p < 0) then p:=-p; q:=-q^ 

13. z;=x; 

14. If ( p 1 |q|to,l) then x:=x + sign(y-x)tol 

else if (p < I" (y-x)q) then x:=x + p/q 

el se xT^ + (y-x)/2; 

15. If (sign f(x) = sign f(y)) then y:=z; 
16., Go to 4; 

21, line 13, Equation (23) with g = 0,... 

27, line 11-23, 

11. k:=k+l 

12. pr=(z-x)f(x) I , 



CORRECTIONS 



PAGE 2 



< > 

13. q:=f{z) - f(x) . ,^ 

1 

14., If (p t< o) then p:=-p; q:=-q; 



18. 



15; z:=x; ^ 

16. m:=m+l 

17. 11^-14) then if (16|y-x| > length) 



n 



■J 



then x:=x+g<y-x); 30 to 20; 
else m:=0; length: = |y-x[ - 



If ( p <^ |q|tol) then x:=x + sign(y-x)tol 
else if. (p < |^y-x)q) 

then x:=x + p/q 

/ • 1 ' 

^' else x:=x + g^y-x) 

19. If_ (sign f(x) =. sign f(y)) then y:=z; 

20. Go to 6 

f 
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